To Do:
library(ggplot2) #for plotting
library(tidyverse)
## ── Attaching packages ────────────────────────────────────────────────────────────────────────────────────────── tidyverse 1.3.0 ──
## ✓ tibble 3.0.1 ✓ dplyr 1.0.0
## ✓ tidyr 1.0.2 ✓ stringr 1.4.0
## ✓ readr 1.3.1 ✓ forcats 0.5.0
## ✓ purrr 0.3.4
## ── Conflicts ───────────────────────────────────────────────────────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(dplyr)
library(GGally)
## Registered S3 method overwritten by 'GGally':
## method from
## +.gg ggplot2
library(gridExtra)
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
The Pinery dataset was compiled in excel from monthly mean temperature and total precipitation data using open access Environment Canada records for stations within 30km of the Pinery with at least 10 consecutive years of data (http://climate.weather.gc.ca/historical_data/search_historic_data_e.html).
The London datasets are homogenized data from Environment Canada that were converted to tidy format in excel.
pinery <- read.csv("data/pinery_climate.csv")
pinery$station <- as.factor(pinery$station) #Make a factor
str(pinery)
## 'data.frame': 2931 obs. of 5 variables:
## $ year : int 1882 1882 1882 1882 1882 1882 1883 1883 1883 1883 ...
## $ month : int 7 8 9 10 11 12 4 5 6 7 ...
## $ mean_temp : num NA NA 16.1 11.4 2.9 -3.9 NA NA NA NA ...
## $ total_precip: num 15.5 182.6 57.9 41.7 28.7 ...
## $ station : Factor w/ 9 levels "Arkona","Centralia A",..: 1 1 1 1 1 1 8 8 8 8 ...
london_precip <- read.csv("data/london_precip.csv")
london_temp <- read.csv("data/london_temp.csv")
london_precip$month <- as.factor(london_precip$month) #Make a factor
london_temp$month <- as.factor(london_temp$month) #Make a factor
climateNA <- read_csv("data/climateNA_monthly.csv")
## Warning: Missing column names filled in: 'X1' [1]
## Parsed with column specification:
## cols(
## X1 = col_double(),
## Year = col_double(),
## Month = col_double(),
## Tmax = col_double(),
## Tmin = col_double(),
## Tave = col_double(),
## PPT = col_double(),
## Rad = col_double(),
## Eref = col_double(),
## CMD = col_double(),
## RH = col_double()
## )
climateNA.temp <- climateNA %>% select(Year, Month, Tave)
climateNA.precip <- climateNA %>% select(Year, Month, PPT)
Here I go through each station for the temp and precip data, removing years where all 12 months aren’t present so that I can later calculate the mean values for each year. I did this manually a while ago and there is probably a more efficient way to do it, however it works and it’s already done so I left it as is.
#Arkona
Arkona <- pinery %>% filter(station == "Arkona") #Filter out station of interest
Arkona_T <- Arkona %>% select(-total_precip) %>% filter(!is.na(Arkona$mean_temp)==T) #remove precipitation data and "NA"
table(Arkona_T$year) #Generate table after cleaning
##
## 1882 1883 1884 1885 1886 1887 1888 1889 1890 1891 1892 1893 1894 1895 1896 1897
## 4 11 11 12 12 12 12 12 12 12 12 5 12 12 12 12
## 1898 1899 1900 1901 1902 1903 1904 1905 1906 1907 1908 1909 1910 1911 1912 1913
## 12 12 10 12 4 12 12 12 12 12 12 12 12 12 12 12
## 1914 1915
## 12 3
Arkona_T2 <- Arkona_T %>% filter(year != 1882 & year != 1883 &
year != 1884 & year != 1893 & year != 1900 & year != 1902 & year != 1915) #remove years with < 12 months
table(Arkona_T2$year) #Generate table after cleaning
##
## 1885 1886 1887 1888 1889 1890 1891 1892 1894 1895 1896 1897 1898 1899 1901 1903
## 12 12 12 12 12 12 12 12 12 12 12 12 12 12 12 12
## 1904 1905 1906 1907 1908 1909 1910 1911 1912 1913 1914
## 12 12 12 12 12 12 12 12 12 12 12
Arkona_P <- Arkona %>% select(-mean_temp) %>% filter(!is.na(Arkona$total_precip)==T) #remove temperature data and "NA"
#table(Arkona_P$year) #show how many observations in each year
Arkona_P2 <- Arkona_P %>% filter(year != 1882 & year != 1893 & year != 1902 & year != 1915)
#table(Arkona_P2$year) #Generate table after cleaning
#Centralia A
CentraliaA <- pinery %>% filter(station == "Centralia A")
CentraliaA_T <- CentraliaA %>% select(-total_precip) %>% filter(!is.na(CentraliaA$mean_temp)==T) #remove "NA"
#table(CentraliaA_T$year)
CentraliaA_T2 <- CentraliaA_T %>% filter(year != 1945 & year != 1946 & year != 1947 & year != 1966)
#table(CentraliaA_T2$year)
CentraliaA_P <- CentraliaA %>% select(-mean_temp) %>% filter(!is.na(CentraliaA$total_precip)==T) #remove "NA"
#table(CentraliaA_P$year)
CentraliaA_P2 <- CentraliaA_P %>% filter(year != 1946 & year != 1947 & year != 1966)
#table(CentraliaA_P2$year)
#Exeter
Exeter <- pinery %>% filter(station == "Exeter")
Exeter_T <- Exeter %>% select(-total_precip) %>% filter(!is.na(Exeter$mean_temp)==T) #remove "NA"
#table(Exeter_T$year)
Exeter_T2 <- Exeter_T %>% filter(year != 1967 & year != 1976 & year != 1985)
#table(Exeter_T2$year)
Exeter_P <- Exeter %>% select(-mean_temp) %>% filter(!is.na(Exeter$total_precip)==T) #remove "NA"
#table(Exeter_P$year)
Exeter_P2 <- Exeter_P %>% filter(year != 1961 & year != 1962 & year != 1964 & year != 1976 & year != 1985)
#table(Exeter_P2$year)
#Forest
Forest <- pinery %>% filter(station == "Forest")
Forest_T <- Forest %>% select(-total_precip) %>% filter(!is.na(Forest$mean_temp)==T) #remove "NA"
#table(Forest_T$year)
Forest_T2 <- Forest_T %>% filter(year != 1924 & year != 1940 & year != 1962 & year != 1963)
#table(Forest_T2$year)
Forest_P <- Forest %>% select(-mean_temp) %>% filter(!is.na(Forest$total_precip)==T) #remove "NA"
#table(Forest_P$year)
Forest_P2 <- Forest_P %>% filter(year > 1933 & year != 1937 & year!= 1943 & year != 1955 & year != 1960 & year < 1960)
#table(Forest_P2$year)
#Lucan
Lucan <- pinery %>% filter(station == "Lucan")
Lucan_T <- Lucan %>% select(-total_precip) %>% filter(!is.na(Lucan$mean_temp)==T) #remove "NA"
#table(Lucan_T$year)
Lucan_T2 <- Lucan_T %>% filter(year != 1915 & year != 1924 & year != 1948 & year != 1949 & year != 1950
& year != 1951 & year != 1952 & year != 1953 & year != 1954 & year != 1955
& year != 1956 & year != 1957 & year != 1958 & year != 1959 & year != 1961)
#table(Lucan_T2$year)
Lucan_P <- Lucan %>% select(-mean_temp) %>% filter(!is.na(Lucan$total_precip)==T) #remove "NA"
#table(Lucan_P$year)
Lucan_P2 <- Lucan_P %>% filter(year != 1915 & year != 1923 & year != 1943 & year < 1948)
#table(Lucan_P2$year)
#Pinery
Pinery <- pinery %>% filter(station == "Pinery")
Pinery_T <- Pinery %>% select(-total_precip) %>% filter(!is.na(Pinery$mean_temp)==T) #remove "NA"
#table(Pinery_T$year)
Pinery_T2 <- Pinery_T %>% filter(year != 1979 & year != 1984)
#table(Pinery_T2$year)
Pinery_P <- Pinery %>% select(-mean_temp) %>% filter(!is.na(Pinery$total_precip)==T) #remove "NA"
#table(Pinery_P$year)
Pinery_P2 <- Pinery_P %>% filter(year != 1979 & year != 1984)
#table(Pinery_P2$year)
#Strathrow-Mullifary
Strathroy <- pinery %>% filter(station == "Strathroy-Mullifary")
Strathroy_T <- Strathroy %>% select(-total_precip) %>% filter(!is.na(Strathroy$mean_temp)==T) #remove "NA"
#table(Strathroy_T$year)
Strathroy_T2 <- Strathroy_T %>% filter(year != 2009 & year != 2010 & year != 2013 & year != 2014 & year != 2015 & year < 2017)
#table(Strathroy_T2$year)
Strathroy_P <- Strathroy %>% select(-mean_temp) %>% filter(!is.na(Strathroy$total_precip)==T) #remove "NA"
#table(Strathroy_P$year)
Strathroy_P2 <- Strathroy_P %>% filter(year != 2009 & year != 2010 & year != 2012 & year != 2013
& year != 2014 & year != 2015 & year < 2017)
#table(Strathroy_P2$year)
#Thedford (no temp)
Thedford <- pinery %>% filter(station == "Thedford")
Thedford_P <- Thedford %>% select(-mean_temp) %>% filter(!is.na(Thedford$total_precip)==T) #remove "NA"
#table(Thedford_P$year)
Thedford_P2 <- Thedford_P %>% filter(year != 1883 & year != 1897)
#table(Thedford_P2$year)
#Watford (can't use precipitation data)
Watford <- pinery %>% filter(station == "Watford")
Watford_T <- Watford %>% select(-total_precip) %>% filter(!is.na(Watford$mean_temp)==T) #remove "NA"
#table(Watford_T$year)
Watford_T2 <- Watford_T %>% filter(year == 1927 | year == 1928)
#table(Watford_T2$year)
Watford_P <- Watford %>% select(-mean_temp) %>% filter(!is.na(Watford$total_precip)==T) #remove "NA"
#table(Watford_P$year)
Re-combine Data
pinery_temp <- rbind(Arkona_T2, CentraliaA_T2, Exeter_T2, Forest_T2, Lucan_T2, Pinery_T2, Strathroy_T2, Watford_T2)
pinery_precip <- rbind(Arkona_P2, CentraliaA_P2, Exeter_P2, Forest_P2, Lucan_P2, Pinery_P2, Strathroy_P2, Thedford_P2)
Ensure the cleaning worked by plotting histograms of the number of observations (months) in each year (should only be multiples of 12)
#histogram (temperature)
graph <- ggplot(pinery_temp, aes(x = year)) +
geom_histogram(position = 'identity', binwidth = 1, color = "black", fill = "grey") +
geom_hline(yintercept = 12, linetype="solid", color = "red") + #shows 12 months
geom_hline(yintercept = 24, linetype="solid", color = "red") +
geom_hline(yintercept = 36, linetype="solid", color = "red") +
theme_bw()
graph
#histogram (precipitation)
graph <- ggplot(pinery_precip, aes(x = year)) +
geom_histogram(position = 'identity', binwidth = 1, color = "black", fill = "grey") +
geom_hline(yintercept = 12, linetype="solid", color = "red") + #shows 12 months
geom_hline(yintercept = 24, linetype="solid", color = "red") +
theme_bw()
graph
Calculate one mean temperature and total precipitation for each year/month where more than one station of data is available
#Temperature
str(pinery_temp)
## 'data.frame': 1944 obs. of 4 variables:
## $ year : int 1885 1885 1885 1885 1885 1885 1885 1885 1885 1885 ...
## $ month : int 1 2 3 4 5 6 7 8 9 10 ...
## $ mean_temp: num -9.4 -13.5 -8.9 4.2 12.1 16.4 21.1 17 14.8 8 ...
## $ station : Factor w/ 9 levels "Arkona","Centralia A",..: 1 1 1 1 1 1 1 1 1 1 ...
pinery_temp_mean <- pinery_temp %>%
group_by(year,month) %>%
summarise(mean_temp_year = mean(mean_temp), sd_temp = sd(mean_temp), n = n()) %>%
mutate(se_temp = sd_temp / sqrt(n), lower = (mean_temp_year - se_temp), upper = (mean_temp_year + se_temp)) %>% round(2)
## `summarise()` regrouping output by 'year' (override with `.groups` argument)
pinery_temp_mean$month <- as.factor(pinery_temp_mean$month)
#Precipitation
pinery_precip_mean <- pinery_precip %>%
group_by(year,month) %>%
summarise(total_precip_year = mean(total_precip), sd_precip = sd(total_precip), n = n()) %>%
mutate(se_precip = sd_precip / sqrt(n), lower = total_precip_year - sd_precip, upper = total_precip_year + sd_precip) %>% round(2)
## `summarise()` regrouping output by 'year' (override with `.groups` argument)
head(pinery_precip_mean)
## # A tibble: 6 x 8
## # Groups: year [1]
## year month total_precip_year sd_precip n se_precip lower upper
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1883 1 30.2 NA 1 NA NA NA
## 2 1883 2 94.7 NA 1 NA NA NA
## 3 1883 3 70.6 NA 1 NA NA NA
## 4 1883 4 33 NA 1 NA NA NA
## 5 1883 5 112. NA 1 NA NA NA
## 6 1883 6 158. NA 1 NA NA NA
ggplot(pinery_precip_mean, aes(x = year)) +
geom_histogram(position = 'identity', binwidth = 1, color = "black", fill = "grey")
ggplot(pinery_temp_mean, aes(x = year)) +
geom_histogram(position = 'identity', binwidth = 1, color = "black", fill = "grey")
Fill empty months-years with data from London
p.t <- pinery_temp_mean %>% select(year, month, mean_temp_year)
l.t <- london_temp %>% filter(month == "Jan" |
month == "Feb" |
month == "Mar" |
month == "Apr" |
month == "May" |
month == "Jun" |
month == "Jul" |
month == "Aug" |
month == "Sep" |
month == "Oct" |
month == "Nov" |
month == "Dec") %>%
mutate(month = recode(month, "Jan"="1",
"Feb"="2",
"Mar"="3",
"Apr"="4",
"May"="5",
"Jun"="6",
"Jul"="7",
"Aug"="8",
"Sep"="9",
"Oct"="10",
"Nov"="11",
"Dec"="12"))
PL.temp <- full_join(l.t, p.t) %>%
rename("pinery" = "mean_temp_year") %>%
rename("london" = "mean_temp")
## Joining, by = c("year", "month")
head(PL.temp)
## year month london pinery
## 1 1883 1 NA NA
## 2 1884 1 -10.1 NA
## 3 1885 1 -9.2 -9.4
## 4 1886 1 -7.2 -7.2
## 5 1887 1 -8.3 -8.4
## 6 1888 1 -10.1 -9.6
Select all rows with regional data - source = Pinery Select all rows with NA - source = London
PL.1 <- PL.temp %>% #use london value where pinery is absent
filter(is.na(pinery)) %>%
mutate(temp = london)
PL.1$source <- "London"
PL.2 <- PL.temp %>% # use pinery value where pinery is present
filter(!is.na(pinery)) %>%
mutate(temp = pinery)
PL.2$source <- "Pinery"
PL.temp2 <- rbind(PL.1, PL.2) %>%
arrange(year)
(nrow(PL.2) + nrow(PL.1)) == nrow(PL.temp) #check equal observations
## [1] TRUE
PL.temp2 %>% filter(is.na(temp))
## year month london pinery temp source
## 1 1883 1 NA NA NA London
## 2 1883 2 NA NA NA London
## 3 2014 9 NA NA NA London
Fill leftover NAs with climateNA values
head(climateNA.temp)
## # A tibble: 6 x 3
## Year Month Tave
## <dbl> <dbl> <dbl>
## 1 1901 1 -4.2
## 2 1902 1 -5.6
## 3 1903 1 -5.1
## 4 1904 1 -9.5
## 5 1905 1 -7.5
## 6 1906 1 -0.7
climateNA.temp <- climateNA.temp %>%
rename("month" = "Month") %>%
rename("year" = "Year") %>%
rename("climateNA" = "Tave")
climateNA.temp$month <- as.factor(climateNA.temp$month)
PL.temp3 <- full_join(PL.temp2, climateNA.temp)
## Joining, by = c("year", "month")
PL.temp3 <- PL.temp3[, c(1:4,7,5,6)] #re-order columns
head(PL.temp3)
## year month london pinery climateNA temp source
## 1 1883 1 NA NA NA NA London
## 2 1883 2 NA NA NA NA London
## 3 1883 3 -6.5 NA NA -6.5 London
## 4 1883 4 5.1 NA NA 5.1 London
## 5 1883 5 10.0 NA NA 10.0 London
## 6 1883 6 17.3 NA NA 17.3 London
PL.3 <- PL.temp3 %>%
filter(is.na(temp)) %>%
mutate(temp = climateNA)
PL.3$source <- "ClimateNA"
PL.4 <- PL.temp3 %>%
filter(!is.na(temp))
PL.temp4 <- rbind(PL.3, PL.4) %>%
arrange(year) %>%
filter(year > 1890) #missing values from 1890s remain
(nrow(PL.3) + nrow(PL.4)) == nrow(PL.temp3) #check equal observations
## [1] TRUE
PL.temp4 %>% filter(is.na(temp)) #values still missing
## [1] year month london pinery climateNA temp source
## <0 rows> (or 0-length row.names)
ggplot(PL.temp4, aes(x = year)) +
geom_histogram(position = 'identity', binwidth = 1, color = "black", fill = "grey") +
geom_hline(yintercept = 12, linetype="solid", color = "red") + #shows 12 months
theme_bw()
p.p <- pinery_precip_mean %>% select(year, month, total_precip_year)
l.p <- london_precip %>% filter(month == "Jan" |
month == "Feb" |
month == "Mar" |
month == "Apr" |
month == "May" |
month == "Jun" |
month == "Jul" |
month == "Aug" |
month == "Sep" |
month == "Oct" |
month == "Nov" |
month == "Dec") %>%
mutate(month = recode(month, "Jan"="1",
"Feb"="2",
"Mar"="3",
"Apr"="4",
"May"="5",
"Jun"="6",
"Jul"="7",
"Aug"="8",
"Sep"="9",
"Oct"="10",
"Nov"="11",
"Dec"="12"))
p.p$month <- as.factor(p.p$month)
PL.precip <- full_join(l.p, p.p) %>%
rename("pinery" = "total_precip_year") %>%
rename("london" = "total_precip")
## Joining, by = c("year", "month")
PL.5 <- PL.precip %>% #use london value where pinery is absent
filter(is.na(pinery)) %>%
mutate(precip = london)
PL.5$source <- "London"
PL.6 <- PL.precip %>% # use pinery value where pinery is present
filter(!is.na(pinery)) %>%
mutate(precip = pinery)
PL.6$source <- "Pinery"
PL.precip2 <- rbind(PL.5, PL.6) %>%
arrange(year)
(nrow(PL.5) + nrow(PL.6)) == nrow(PL.precip2) #check equal observations
## [1] TRUE
PL.precip2 %>% filter(is.na(precip)) #values still missing
## year month london pinery precip source
## 1 2009 4 NA NA NA London
## 2 2009 5 NA NA NA London
## 3 2009 6 NA NA NA London
## 4 2009 7 NA NA NA London
## 5 2009 8 NA NA NA London
## 6 2009 9 NA NA NA London
## 7 2009 10 NA NA NA London
## 8 2010 4 NA NA NA London
## 9 2010 5 NA NA NA London
## 10 2010 6 NA NA NA London
## 11 2010 7 NA NA NA London
## 12 2010 8 NA NA NA London
## 13 2010 9 NA NA NA London
## 14 2010 10 NA NA NA London
## 15 2012 4 NA NA NA London
## 16 2012 5 NA NA NA London
## 17 2012 6 NA NA NA London
## 18 2012 7 NA NA NA London
## 19 2012 8 NA NA NA London
## 20 2012 9 NA NA NA London
## 21 2012 10 NA NA NA London
## 22 2013 4 NA NA NA London
## 23 2013 5 NA NA NA London
## 24 2013 6 NA NA NA London
## 25 2013 7 NA NA NA London
## 26 2013 8 NA NA NA London
## 27 2013 9 NA NA NA London
## 28 2013 10 NA NA NA London
## 29 2014 4 NA NA NA London
## 30 2014 5 NA NA NA London
## 31 2014 6 NA NA NA London
## 32 2014 7 NA NA NA London
## 33 2014 8 NA NA NA London
## 34 2014 9 NA NA NA London
## 35 2014 10 NA NA NA London
## 36 2015 4 NA NA NA London
## 37 2015 5 NA NA NA London
## 38 2015 6 NA NA NA London
## 39 2015 7 NA NA NA London
## 40 2015 8 NA NA NA London
## 41 2015 9 NA NA NA London
## 42 2015 10 NA NA NA London
## 43 2017 4 NA NA NA London
## 44 2017 5 NA NA NA London
## 45 2017 6 NA NA NA London
## 46 2017 7 NA NA NA London
## 47 2017 8 NA NA NA London
## 48 2017 9 NA NA NA London
## 49 2017 10 NA NA NA London
## 50 2017 11 NA NA NA London
## 51 2017 12 NA NA NA London
Fill leftover NAs with climateNA values
head(climateNA.temp)
## # A tibble: 6 x 3
## year month climateNA
## <dbl> <fct> <dbl>
## 1 1901 1 -4.2
## 2 1902 1 -5.6
## 3 1903 1 -5.1
## 4 1904 1 -9.5
## 5 1905 1 -7.5
## 6 1906 1 -0.7
climateNA.precip <- climateNA.precip %>%
rename("month" = "Month") %>%
rename("year" = "Year") %>%
rename("climateNA" = "PPT")
climateNA.precip$month <- as.factor(climateNA.precip$month)
PL.precip3 <- full_join(PL.precip2, climateNA.precip)
## Joining, by = c("year", "month")
PL.precip3 <- PL.precip3[, c(1:4,7,5,6)] #re-order columns
head(PL.temp3)
## year month london pinery climateNA temp source
## 1 1883 1 NA NA NA NA London
## 2 1883 2 NA NA NA NA London
## 3 1883 3 -6.5 NA NA -6.5 London
## 4 1883 4 5.1 NA NA 5.1 London
## 5 1883 5 10.0 NA NA 10.0 London
## 6 1883 6 17.3 NA NA 17.3 London
PL.7 <- PL.precip3 %>%
filter(is.na(precip)) %>%
mutate(precip = climateNA)
PL.7$source <- "ClimateNA"
PL.8 <- PL.precip3 %>%
filter(!is.na(precip))
PL.7
## year month london pinery climateNA precip source
## 1 2009 4 NA NA 132 132 ClimateNA
## 2 2009 5 NA NA 67 67 ClimateNA
## 3 2009 6 NA NA 108 108 ClimateNA
## 4 2009 7 NA NA 72 72 ClimateNA
## 5 2009 8 NA NA 90 90 ClimateNA
## 6 2009 9 NA NA 40 40 ClimateNA
## 7 2009 10 NA NA 116 116 ClimateNA
## 8 2010 4 NA NA 64 64 ClimateNA
## 9 2010 5 NA NA 103 103 ClimateNA
## 10 2010 6 NA NA 132 132 ClimateNA
## 11 2010 7 NA NA 100 100 ClimateNA
## 12 2010 8 NA NA 36 36 ClimateNA
## 13 2010 9 NA NA 76 76 ClimateNA
## 14 2010 10 NA NA 57 57 ClimateNA
## 15 2012 4 NA NA 58 58 ClimateNA
## 16 2012 5 NA NA 50 50 ClimateNA
## 17 2012 6 NA NA 52 52 ClimateNA
## 18 2012 7 NA NA 67 67 ClimateNA
## 19 2012 8 NA NA 70 70 ClimateNA
## 20 2012 9 NA NA 73 73 ClimateNA
## 21 2012 10 NA NA 128 128 ClimateNA
## 22 2013 4 NA NA 115 115 ClimateNA
## 23 2013 5 NA NA 71 71 ClimateNA
## 24 2013 6 NA NA 131 131 ClimateNA
## 25 2013 7 NA NA 91 91 ClimateNA
## 26 2013 8 NA NA 85 85 ClimateNA
## 27 2013 9 NA NA 51 51 ClimateNA
## 28 2013 10 NA NA 123 123 ClimateNA
## 29 2014 4 NA NA 81 81 ClimateNA
## 30 2014 5 NA NA 81 81 ClimateNA
## 31 2014 6 NA NA 88 88 ClimateNA
## 32 2014 7 NA NA 90 90 ClimateNA
## 33 2014 8 NA NA 103 103 ClimateNA
## 34 2014 9 NA NA 95 95 ClimateNA
## 35 2014 10 NA NA 84 84 ClimateNA
## 36 2015 4 NA NA 74 74 ClimateNA
## 37 2015 5 NA NA 95 95 ClimateNA
## 38 2015 6 NA NA 118 118 ClimateNA
## 39 2015 7 NA NA 65 65 ClimateNA
## 40 2015 8 NA NA 79 79 ClimateNA
## 41 2015 9 NA NA 59 59 ClimateNA
## 42 2015 10 NA NA 68 68 ClimateNA
## 43 2017 4 NA NA 118 118 ClimateNA
## 44 2017 5 NA NA 122 122 ClimateNA
## 45 2017 6 NA NA 116 116 ClimateNA
## 46 2017 7 NA NA 75 75 ClimateNA
## 47 2017 8 NA NA 82 82 ClimateNA
## 48 2017 9 NA NA 86 86 ClimateNA
## 49 2017 10 NA NA 107 107 ClimateNA
## 50 2017 11 NA NA 110 110 ClimateNA
## 51 2017 12 NA NA 50 50 ClimateNA
## 52 2018 1 NA NA 66 66 ClimateNA
## 53 2018 2 NA NA 107 107 ClimateNA
## 54 2018 3 NA NA 68 68 ClimateNA
## 55 2018 4 NA NA 102 102 ClimateNA
## 56 2018 5 NA NA 67 67 ClimateNA
## 57 2018 6 NA NA 97 97 ClimateNA
## 58 2018 7 NA NA 89 89 ClimateNA
## 59 2018 8 NA NA 180 180 ClimateNA
## 60 2018 9 NA NA 71 71 ClimateNA
## 61 2018 10 NA NA 100 100 ClimateNA
## 62 2018 11 NA NA 136 136 ClimateNA
## 63 2018 12 NA NA 58 58 ClimateNA
PL.precip4 <- rbind(PL.7, PL.8) %>%
arrange(year) %>%
filter(year >=1890)
(nrow(PL.7) + nrow(PL.8)) == nrow(PL.precip3) #check equal observations
## [1] TRUE
PL.precip4 %>% filter(is.na(precip)) #values still missing
## [1] year month london pinery climateNA precip source
## <0 rows> (or 0-length row.names)
ggplot(PL.precip4, aes(x = year)) +
geom_histogram(position = 'identity', binwidth = 1, color = "black", fill = "grey") +
geom_hline(yintercept = 12, linetype="solid", color = "red") + #shows 12 months
theme_bw()
study_t <- pinery_temp_mean %>%
group_by(month) %>%
summarize(mean_temp = mean(mean_temp_year) %>% round(2))
## `summarise()` ungrouping output (override with `.groups` argument)
study_t
## # A tibble: 12 x 2
## month mean_temp
## <fct> <dbl>
## 1 1 -5.64
## 2 2 -5.78
## 3 3 -0.64
## 4 4 6.37
## 5 5 12.8
## 6 6 18.2
## 7 7 20.8
## 8 8 19.9
## 9 9 16.2
## 10 10 9.83
## 11 11 3.28
## 12 12 -2.8
study_p <- pinery_precip_mean %>%
group_by(month) %>%
summarize(total_precip= mean(total_precip_year), sd_precip = sd(total_precip_year), n = n()) %>%
mutate(se_precip = sd_precip / sqrt(n), lower = total_precip - se_precip, upper = total_precip + se_precip %>% round(2))
## `summarise()` ungrouping output (override with `.groups` argument)
study_p
## # A tibble: 12 x 7
## month total_precip sd_precip n se_precip lower upper
## <dbl> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 1 80.5 32.6 122 2.95 77.6 83.5
## 2 2 65.4 26.8 122 2.42 63.0 67.8
## 3 3 64.4 26.0 122 2.35 62.1 66.8
## 4 4 72.5 32.0 122 2.90 69.7 75.4
## 5 5 82.0 38.0 122 3.44 78.6 85.5
## 6 6 78.5 36.3 122 3.28 75.2 81.8
## 7 7 81.4 48.6 122 4.40 77.0 85.8
## 8 8 79.3 43.4 122 3.93 75.4 83.3
## 9 9 86.2 51.0 122 4.61 81.6 90.8
## 10 10 82.9 41.0 122 3.72 79.2 86.6
## 11 11 89.1 32.2 122 2.91 86.2 92.0
## 12 12 87.2 31.9 122 2.89 84.3 90.1
str(study_t)
## tibble [12 × 2] (S3: tbl_df/tbl/data.frame)
## $ month : Factor w/ 12 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ mean_temp: num [1:12] -5.64 -5.78 -0.64 6.37 12.83 ...
ggplot() +
geom_col(data = study_t, aes(x = month, y = mean_temp), fill = "red3", alpha = 0.7) +
geom_line(data = study_p, aes(x = month, y = total_precip-80), col= "steelblue4", size = 1.5) +
scale_y_continuous("Temperature (°C)", sec.axis = sec_axis(~ . + 80, name = "Precipitation (mm)")) +
geom_ribbon(data = study_p, aes(x = month, ymin = (lower-80), ymax = (upper-80)), alpha = 0.3) +
scale_x_discrete("Month") +
theme_bw() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
p.climograph <- ggplot() +
geom_boxplot(data = pinery_temp_mean, aes(x = month, y = mean_temp_year, group = month), fill = "red3", alpha = 0.7) +
geom_line(data = study_p, aes(x = month, y = total_precip-75), col= "steelblue4", size = 1) +
geom_ribbon(data = study_p, aes(x = month, ymin = (lower-75), ymax = (upper-75)), fill = "steelblue4", alpha = 0.2) +
scale_y_continuous("Temperature (°C)", sec.axis = sec_axis(~ . + 75, name = "Precipitation (mm)")) +
scale_x_discrete("Month") +
theme_bw() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
p.climograph
#ggsave('p.climograph.pdf', p.climograph, units = 'cm', width = 20, height = 12)
london_temp_annual <- london_temp %>% filter(month == "Annual")
london_precip_annual <- london_precip %>% filter(month == "Annual")
pinery_temp_annual <- pinery_temp_mean %>% group_by(year) %>% summarise(mean_temp = mean(mean_temp_year))
## `summarise()` ungrouping output (override with `.groups` argument)
pinery_precip_annual <- pinery_precip_mean %>% group_by(year) %>% summarise(total_precip = sum(total_precip_year))
## `summarise()` ungrouping output (override with `.groups` argument)
london_temp_2 <- london_temp_annual %>% select(-month)
london_temp_2$location <- "London"
pinery_temp_annual$location <- "Pinery"
annual_temp <- rbind(london_temp_2,pinery_temp_annual) #join data sets
london_precip_2 <- london_precip_annual %>% select(-month)
london_precip_2$location <- "London"
pinery_precip_annual$location <- "Pinery"
annual_precip <- rbind(london_precip_2,pinery_precip_annual) #join data sets
Plot Temperature
annual_temp_plot <- ggplot() +
# London Data
geom_line(data=london_temp_annual, aes(x = year, y = mean_temp), colour = "lightseagreen", linetype = "dashed", alpha = 0.8) +
# Pinery Data
geom_line(data=pinery_temp_annual, aes(x = year, y = mean_temp), colour="lightsalmon4") +
#Aesthetics
xlab("Year") +
ylab("Mean Annual Temperature (°C)") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
#ggtitle("Average Annual Temperature in \n Southwestern Ontario from 1883 to 2017") +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
annual_temp_plot
## Warning: Removed 1 row(s) containing missing values (geom_path).
#ggsave('annual_temp.pdf', annual_temp_plot, units = 'cm', width = 15, height = 10)
Plot recent temperature
pinery_temp_annual_recent <- pinery_temp_annual %>% filter(year > 1949)
annual_temp_plot2 <- ggplot() +
# Pinery Data
geom_line(data=pinery_temp_annual_recent, aes(x = year, y = mean_temp)) +
geom_point(data=pinery_temp_annual_recent, aes(x = year, y = mean_temp), size = 1) +
#Aesthetics
xlab("Year") +
ylab("Temperature (°C)") +
scale_x_continuous(breaks=seq(1950,2020,10)) +
#ggtitle("Average Annual Temperature in \n Southwestern Ontario from 1883 to 2017") +
scale_colour_manual(name="Location", breaks = c("London", "Pinery"),
values = c("blue","darkgreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
#ggsave('annual_temp.pdf', annual_temp_plot, units = 'cm', width = 20, height = 12)
annual_temp_plot2
Plot Precipitation
annual_precip_plot <- ggplot() +
# London Data
geom_line(data=london_precip_annual, aes(x = year, y = total_precip), colour = "lightseagreen", linetype = "dashed", alpha = 0.8) +
# geom_smooth(data=london_precip_annual, aes(x = year, y = total_precip), se = TRUE,
# method = 'lm', size = 0.8, colour = "black", fill = "blue") +
# Pinery Data
geom_line(data=pinery_precip_annual, aes(x = year, y = total_precip, group = 1), colour="lightsalmon4")+
# geom_smooth(data=pinery_precip_annual, aes(x = year, y = total_precip), se = TRUE,
# method = 'lm', size = 0.8, colour = "black", fill = "darkgreen") +
#Aesthetics
xlab("Year") +
ylab("Total Annual Precipitation (mm)") +
scale_x_continuous(breaks=seq(1880,2020,20)) +
#ggtitle("Annual Precipitation in Southwestern Ontario \nfrom 1883 to 2017") +
scale_colour_manual(name="Location", breaks = c("London", "Pinery"),
values = c("blue","darkgreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
annual_precip_plot
## Warning: Removed 17 row(s) containing missing values (geom_path).
#ggsave('annual_precip.pdf', annual_precip_plot, units = 'cm', width = 15, height = 10)
annual_plot <- grid.arrange(annual_temp_plot, annual_precip_plot, ncol= 2)
## Warning: Removed 1 row(s) containing missing values (geom_path).
## Warning: Removed 17 row(s) containing missing values (geom_path).
#ggsave('annual_plot.pdf', annual_plot, units = 'cm', width = 23, height = 8)
winter_t <- pinery_temp %>% filter(month == 1 | month == 2 | month == 3)
winter_t$season <- "winter"
spring_t <- pinery_temp %>% filter(month == 4 | month == 5 | month == 6)
spring_t$season <- "spring"
summer_t <- pinery_temp %>% filter(month == 7 | month == 8 | month == 9)
summer_t$season <- "summer"
fall_t <- pinery_temp %>% filter(month == 10 | month == 11 | month == 12)
fall_t$season <- "fall"
pinery_temp <- rbind(winter_t, spring_t, summer_t, fall_t)
pinery_temp$season <- as.factor(pinery_temp$season) #Make a factor
str(pinery_temp)
## 'data.frame': 1944 obs. of 5 variables:
## $ year : int 1885 1885 1885 1886 1886 1886 1887 1887 1887 1888 ...
## $ month : int 1 2 3 1 2 3 1 2 3 1 ...
## $ mean_temp: num -9.4 -13.5 -8.9 -7.2 -6.9 -1 -8.4 -5.4 -3.8 -9.6 ...
## $ station : Factor w/ 9 levels "Arkona","Centralia A",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ season : Factor w/ 4 levels "fall","spring",..: 4 4 4 4 4 4 4 4 4 4 ...
pinery_temp_station <- pinery_temp %>%
group_by(year, season, station) %>%
summarise(mean_temp = mean(mean_temp))
## `summarise()` regrouping output by 'year', 'season' (override with `.groups` argument)
pinery_temp_season <- pinery_temp_station %>%
group_by(year, season) %>%
summarise(mean_temp_season = mean(mean_temp) %>% round(2), sd_temp = sd(mean_temp) %>% round(2), n = n()) %>%
mutate(se_temp = (sd_temp / sqrt(n)) %>% round(2), lower = mean_temp_season - se_temp, upper = mean_temp_season + se_temp)
## `summarise()` regrouping output by 'year' (override with `.groups` argument)
head(pinery_temp_season)
## # A tibble: 6 x 8
## # Groups: year [2]
## year season mean_temp_season sd_temp n se_temp lower upper
## <int> <fct> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 1885 fall 2.7 NA 1 NA NA NA
## 2 1885 spring 10.9 NA 1 NA NA NA
## 3 1885 summer 17.6 NA 1 NA NA NA
## 4 1885 winter -10.6 NA 1 NA NA NA
## 5 1886 fall 1.33 NA 1 NA NA NA
## 6 1886 spring 12.4 NA 1 NA NA NA
winter_p <- pinery_precip %>% filter(month == 1 | month == 2 | month == 3)
winter_p$season <- "winter"
spring_p <- pinery_precip %>% filter(month == 4 | month == 5 | month == 6)
spring_p$season <- "spring"
summer_p <- pinery_precip %>% filter(month == 7 | month == 8 | month == 9)
summer_p$season <- "summer"
fall_p <- pinery_precip %>% filter(month == 10 | month == 11 | month == 12)
fall_p$season <- "fall"
pinery_precip <- rbind(winter_p, spring_p, summer_p, fall_p)
pinery_precip$season <- as.factor(pinery_precip$season) #Make a factor
pinery_precip_station <- pinery_precip %>%
group_by(year, season, station) %>%
summarise(total_precip = mean(total_precip))
## `summarise()` regrouping output by 'year', 'season' (override with `.groups` argument)
pinery_precip_season <- pinery_precip_station %>%
group_by(year, season) %>%
summarise(total_precip_season = mean(total_precip) %>% round(2), sd_precip = sd(total_precip) %>% round(2), n = n()) %>%
mutate(se_precip = sd_precip / sqrt(n), lower = total_precip_season - sd_precip, upper = total_precip_season + sd_precip %>% round(2))
## `summarise()` regrouping output by 'year' (override with `.groups` argument)
head(pinery_precip_season)
## # A tibble: 6 x 8
## # Groups: year [2]
## year season total_precip_season sd_precip n se_precip lower upper
## <int> <fct> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 1883 fall 69.3 NA 1 NA NA NA
## 2 1883 spring 101. NA 1 NA NA NA
## 3 1883 summer 83.8 NA 1 NA NA NA
## 4 1883 winter 65.2 NA 1 NA NA NA
## 5 1884 fall 91.5 1.01 2 0.714 90.5 92.5
## 6 1884 spring 49.1 1.51 2 1.07 47.6 50.6
winter_tl <- london_temp %>% filter(month == "Jan" | month == "Feb" | month == "Mar")
winter_tl$season <- "winter"
spring_tl <- london_temp %>% filter(month == "Apr" | month == "May" | month == "Jun")
spring_tl$season <- "spring"
summer_tl <- london_temp %>% filter(month == "Jul" | month == "Aug" | month == "Sep")
summer_tl$season <- "summer"
fall_tl <- london_temp %>% filter(month == "Oct" | month == "Nov" | month == "Dec")
fall_tl$season <- "fall"
london_temp <- rbind(winter_tl, spring_tl, summer_tl, fall_tl)
london_temp$season <- as.factor(london_temp$season) #Make a factor
str(london_temp)
## 'data.frame': 1620 obs. of 4 variables:
## $ year : int 1883 1884 1885 1886 1887 1888 1889 1890 1891 1892 ...
## $ month : Factor w/ 17 levels "Annual","Apr",..: 7 7 7 7 7 7 7 7 7 7 ...
## $ mean_temp: num NA -10.1 -9.2 -7.2 -8.3 -10.1 -3 -0.7 -4.3 -7.9 ...
## $ season : Factor w/ 4 levels "fall","spring",..: 4 4 4 4 4 4 4 4 4 4 ...
london_temp_season <- london_temp %>%
group_by(year, season) %>%
summarise(mean_temp_season = mean(mean_temp))
## `summarise()` regrouping output by 'year' (override with `.groups` argument)
winter_pl <- london_precip %>% filter(month == "Jan" | month == "Feb" | month == "Mar")
winter_pl$season <- "winter"
spring_pl <- london_precip %>% filter(month == "Apr" | month == "May" | month == "Jun")
spring_pl$season <- "spring"
summer_pl <- london_precip %>% filter(month == "Jul" | month == "Aug" | month == "Sep")
summer_pl$season <- "summer"
fall_pl <- london_precip %>% filter(month == "Oct" | month == "Nov" | month == "Dec")
fall_pl$season <- "fall"
london_precip <- rbind(winter_pl, spring_pl, summer_pl, fall_pl)
london_precip$season <- as.factor(london_precip$season) #Make a factor
str(london_precip)
## 'data.frame': 1620 obs. of 4 variables:
## $ year : int 1883 1884 1885 1886 1887 1888 1889 1890 1891 1892 ...
## $ month : Factor w/ 17 levels "Annual","Apr",..: 7 7 7 7 7 7 7 7 7 7 ...
## $ total_precip: num NA 106.2 96.1 142.5 86.9 ...
## $ season : Factor w/ 4 levels "fall","spring",..: 4 4 4 4 4 4 4 4 4 4 ...
london_precip_season <- london_precip %>%
group_by(year, season) %>%
summarise(total_precip_season = mean(total_precip))
## `summarise()` regrouping output by 'year' (override with `.groups` argument)
Make a factor
london_temp_season$season <- factor(london_temp_season$season, levels = c("winter", "spring", "summer", "fall"))
pinery_temp_season$season <- factor(pinery_temp_season$season, levels = c("winter", "spring", "summer", "fall"))
pinery_temp_season$location <- "Pinery"
london_temp_season$location <- "London"
season_temp <- rbind(pinery_temp_season, london_temp_season)
pinery_precip_season$location <- "Pinery"
london_precip_season$location <- "London"
season_precip <- rbind(pinery_precip_season, london_precip_season)
season_temp$location <- factor(season_temp$location, levels = c("Pinery", "London"))
season_precip$location <- factor(season_precip$location, levels = c("Pinery", "London"))
Plot Temperature
w_temp <- ggplot() +
geom_line(data=(season_temp %>% filter(season == "winter")), aes(x = year, y = mean_temp_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature (°C)") + ggtitle("Winter") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
sp_temp <- ggplot() +
geom_line(data=(season_temp %>% filter(season == "spring")), aes(x = year, y = mean_temp_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature (°C)") + ggtitle("Spring") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
su_temp <- ggplot() +
geom_line(data=(season_temp %>% filter(season == "summer")), aes(x = year, y = mean_temp_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature (°C)") + ggtitle("Summer") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
fa_temp <- ggplot() +
geom_line(data=(season_temp %>% filter(season == "fall")), aes(x = year, y = mean_temp_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature (°C)") + ggtitle("Fall") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
seasonal_temp_plot <- grid.arrange(w_temp, sp_temp, su_temp, fa_temp, nrow = 2, ncol = 2)
## Warning: Removed 1 row(s) containing missing values (geom_path).
#ggsave('seasonal_temperature.pdf', seasonal_temp_plot, units = 'cm', width = 20, height = 12)
Plot Precipitation
w_precip <- ggplot() +
geom_line(data=(season_precip %>% filter(season == "winter")), aes(x = year, y = total_precip_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Precipitation (mm)") + ggtitle("Winter") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
sp_precip <- ggplot() +
geom_line(data=(season_precip %>% filter(season == "spring")), aes(x = year, y = total_precip_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Precipitation (mm)") + ggtitle("Spring") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
su_precip <- ggplot() +
geom_line(data=(season_precip %>% filter(season == "summer")), aes(x = year, y = total_precip_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Precipitation (mm)") + ggtitle("Summer") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
fa_precip <- ggplot() +
geom_line(data=(season_precip %>% filter(season == "fall")), aes(x = year, y = total_precip_season, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Precipitation (mm)") + ggtitle("Fall") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
scale_color_manual(values = c("lightsalmon4", "lightseagreen")) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
seasonal_precip_plot <- grid.arrange(w_precip, sp_precip, su_precip, fa_precip, nrow = 2, ncol = 2)
## Warning: Removed 1 row(s) containing missing values (geom_path).
## Warning: Removed 15 row(s) containing missing values (geom_path).
## Warning: Removed 15 row(s) containing missing values (geom_path).
## Warning: Removed 16 row(s) containing missing values (geom_path).
#ggsave('seasonal_precipitation.pdf', seasonal_precip_plot, units = 'cm', width = 20, height = 12)
https://tamino.wordpress.com/2018/07/29/why-use-temperature-anomaly/
Average defined as 1951-1980 (NASA standard)
pinery_average <- season_temp %>% filter(location == "Pinery") %>% filter(year > 1950 & year < 1981)
pinery_average <- pinery_average %>% group_by(season) %>% summarise(pinery_average = mean(mean_temp_season))
## `summarise()` ungrouping output (override with `.groups` argument)
london_average <- season_temp %>% filter(location == "London") %>% filter(year > 1950 & year < 1981)
london_average <- london_average %>% group_by(season) %>% summarise(london_average = mean(mean_temp_season))
## `summarise()` ungrouping output (override with `.groups` argument)
pinery_anomaly <- right_join(season_temp, pinery_average)
## Joining, by = "season"
pinery_anomaly <- mutate(pinery_anomaly, anomaly = mean_temp_season - pinery_average) %>% filter(location == "Pinery")
london_anomaly <- right_join(season_temp, london_average)
## Joining, by = "season"
london_anomaly <- mutate(london_anomaly, anomaly = mean_temp_season - london_average) %>% filter(location == "London")
season_anomaly <- full_join(pinery_anomaly, london_anomaly) %>% select(-pinery_average, -london_average)
## Joining, by = c("year", "season", "mean_temp_season", "sd_temp", "n", "se_temp", "lower", "upper", "location", "anomaly")
w_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "winter")), aes(x = year, y = anomaly, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Winter") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
sp_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "spring")), aes(x = year, y = anomaly, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Spring") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
su_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "summer")), aes(x = year, y = anomaly, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Summer") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
fa_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "fall")), aes(x = year, y = anomaly, col = location, linetype = location)) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Fall") +
scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
seasonal_anomaly_plot <- grid.arrange(w_anomaly, sp_anomaly, su_anomaly, fa_anomaly, nrow = 2, ncol = 2)
## Warning: Removed 1 row(s) containing missing values (geom_path).
#ggsave('seasonal_anomaly.pdf', seasonal_anomaly_plot, units = 'cm', width = 20, height = 12)
Plot recent temperature anomaly
w_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "winter", year > 1949, location == "Pinery")), aes(x = year, y = anomaly)) +
geom_point(data=(season_anomaly %>% filter(season == "winter", year > 1949, location == "Pinery")), aes(x = year, y = anomaly), size = 0.5) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Winter") +
scale_x_continuous(limits = c(1950, 2020), breaks=seq(1950,2020,10)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
sp_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "spring", year > 1949, location == "Pinery")), aes(x = year, y = anomaly)) +
geom_point(data=(season_anomaly %>% filter(season == "spring", year > 1949, location == "Pinery")), aes(x = year, y = anomaly), size = 0.5) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Spring") +
scale_x_continuous(limits = c(1950, 2020), breaks=seq(1950,2020,10)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
su_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "summer", year > 1949, location == "Pinery")), aes(x = year, y = anomaly)) +
geom_point(data=(season_anomaly %>% filter(season == "summer", year > 1949, location == "Pinery")), aes(x = year, y = anomaly), size = 0.5) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Summer") +
scale_x_continuous(limits = c(1950, 2020), breaks=seq(1950,2020,10)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
fa_anomaly <- ggplot() +
geom_line(data=(season_anomaly %>% filter(season == "fall", year > 1949, location == "Pinery")), aes(x = year, y = anomaly)) +
geom_point(data=(season_anomaly %>% filter(season == "fall", year > 1949, location == "Pinery")), aes(x = year, y = anomaly), size = 0.5) +
#Aesthetics
xlab("Year") + ylab("Temperature Anomaly (°C)") + ggtitle("Fall") +
scale_x_continuous(limits = c(1950, 2020), breaks=seq(1950,2020,10)) +
theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.position = "none")
seasonal_anomaly_plot <- grid.arrange(w_anomaly, sp_anomaly, su_anomaly, fa_anomaly, nrow = 2, ncol = 2)
#write.csv(pinery_temp_mean, "outputs/monthly_temperature.csv")
#write.csv(pinery_precip_mean, "outputs/monthly_precipitation.csv")
#write.csv(season_temp, file = "outputs/season_temperature.csv")
#write.csv(season_precip, file = "outputs/season_precipitation.csv")
#write.csv(pinery_temp_annual, file = "outputs/annual_temp.csv")
Sys.time()
## [1] "2020-08-05 16:30:07 PDT"
git2r::repository()
## Local: master /Users/JenBaron/Documents/UWO/UWO NSERC/Statistical Analysis/Climate/Climate_Pinery
## Remote: master @ origin (https://github.com/JenBaron/Climate_Pinery.git)
## Head: [a85d3a9] 2020-08-05: Fill missing values with London & ClimateNA
sessionInfo()
## R version 4.0.0 (2020-04-24)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Catalina 10.15.5
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_CA.UTF-8/en_CA.UTF-8/en_CA.UTF-8/C/en_CA.UTF-8/en_CA.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] gridExtra_2.3 GGally_1.5.0 forcats_0.5.0 stringr_1.4.0
## [5] dplyr_1.0.0 purrr_0.3.4 readr_1.3.1 tidyr_1.0.2
## [9] tibble_3.0.1 tidyverse_1.3.0 ggplot2_3.3.1
##
## loaded via a namespace (and not attached):
## [1] tidyselect_1.1.0 xfun_0.14 haven_2.2.0 lattice_0.20-41
## [5] colorspace_1.4-1 vctrs_0.3.1 generics_0.0.2 htmltools_0.4.0
## [9] yaml_2.2.1 utf8_1.1.4 rlang_0.4.6 pillar_1.4.4
## [13] glue_1.4.1 withr_2.2.0 DBI_1.1.0 RColorBrewer_1.1-2
## [17] dbplyr_1.4.3 modelr_0.1.6 readxl_1.3.1 plyr_1.8.6
## [21] lifecycle_0.2.0 munsell_0.5.0 gtable_0.3.0 cellranger_1.1.0
## [25] rvest_0.3.5 evaluate_0.14 labeling_0.3 knitr_1.28
## [29] fansi_0.4.1 broom_0.5.6 Rcpp_1.0.4.6 scales_1.1.1
## [33] backports_1.1.7 jsonlite_1.6.1 farver_2.0.3 fs_1.4.1
## [37] hms_0.5.3 digest_0.6.25 stringi_1.4.6 grid_4.0.0
## [41] cli_2.0.2 tools_4.0.0 magrittr_1.5 crayon_1.3.4
## [45] pkgconfig_2.0.3 ellipsis_0.3.1 xml2_1.3.2 reprex_0.3.0
## [49] lubridate_1.7.8 reshape_0.8.8 assertthat_0.2.1 rmarkdown_2.1
## [53] httr_1.4.1 rstudioapi_0.11 R6_2.4.1 git2r_0.26.1
## [57] nlme_3.1-147 compiler_4.0.0